#### Replication code for "The Black Market Blues: The Political Costs of Illicit Currency Markets"

#### Presidential approval #####

rm(list = ls())
library(dplyr)
library(ggplot2)
library(reshape2)
library(stargazer)
library(readstata13)


rm(list = ls())

setwd('/Users/lschiume/Dropbox/Argentina Black Market/JOP Final Submission')
ts<-read.dta13('Black_Market_Blues_Time-series.dta')
ts$id<-seq(1,dim(ts)[1],1)

ts$cepo<-ifelse(ts$ceposum5>1,1,0)
ts$bmp<-ts$bmp*100
lag<-ts[,c( 'id','approval','bmpgfd','cepo1','rer1', 'gdpgrowth','bppinflation','time',
            'electioncycle','electioncycle2')]


colnames(lag)[2:length(lag)]<-paste('l.',colnames(lag)[2:length(lag)],sep='')
lag$id<-lag$id+1
ts<-merge(ts,lag,by.x='id',by.y='id',all.x=T)
ts$d.approval<-ts$approval-ts$l.approval
ts$d.bmpgfd<-ts$bmpgfd-ts$l.bmpgfd
ts$d.cepo1<-ts$cepo1-ts$l.cepo1
ts$d.bppinflation<-ts$bppinflation-ts$l.bppinflation
ts$d.gdpgrowth<-ts$gdpgrowth-ts$l.gdpgrowth
ts$d.rer1<-ts$rer-ts$l.rer


lag<-ts[,c( 'id','d.approval','d.bmpgfd','d.cepo1','d.rer1', 'd.gdpgrowth','d.bppinflation')]

colnames(lag)[2:length(lag)]<-paste('l.',colnames(lag)[2:length(lag)],sep='')
lag$id<-lag$id+1
ts<-merge(ts,lag,by.x='id',by.y='id',all.x=T)


# Table 1
models<-list()

models$m1<-lm(d.approval~l.approval +
                l.d.bmpgfd +l.bmpgfd + l.d.cepo1+  l.cepo1 +l.d.rer1 + l.rer1 +
                l.d.gdpgrowth + l.gdpgrowth + l.d.bppinflation + l.bppinflation,data=ts)

models$m2<-lm(d.approval~l.approval +
              l.d.bmpgfd +l.bmpgfd + l.d.cepo1+  l.cepo1 +l.d.rer1 + l.rer1 +
                l.d.gdpgrowth + l.gdpgrowth + l.d.bppinflation + l.bppinflation + time,data=ts)


models$m3<-lm(d.approval~l.approval +
                l.d.bmpgfd +l.bmpgfd + l.d.cepo1+  l.cepo1 +l.d.rer1 + l.rer1 +
                l.d.gdpgrowth + l.gdpgrowth + l.d.bppinflation + l.bppinflation + electioncycle+electioncycle2,data=ts)


stargazer(models,type='text',keep.stat=c('n'))


# Figure 1

alpha1<-coef(models$m1)[2]
alpha2<-coef(models$m2)[2]
alpha3<-coef(models$m3)[2]
beta1.bmp1<-coef(models$m1)[4]
beta1.bmp2<-coef(models$m2)[4]
beta1.bmp3<-coef(models$m3)[4]

change.bmp<-10
lrm.bmp1<--(beta1.bmp1/alpha1)*change.bmp
lrm.bmp2<--(beta1.bmp2/alpha2)*change.bmp
lrm.bmp3<--(beta1.bmp3/alpha3)*change.bmp
model1<-matrix(nrow=100,ncol=1,NA)
model2<-matrix(nrow=100,ncol=1,NA)
model3<-matrix(nrow=100,ncol=1,NA)



for(i in 1:100){
  if (i==1){
    bmp1<-beta1.bmp1*change.bmp
    bmp2<-beta1.bmp2*change.bmp
    bmp3<-beta1.bmp3*change.bmp
  }else{
    bmp1<-bmp1*(1+alpha1)
    bmp2<-bmp2*(1+alpha1)
    bmp3<-bmp3*(1+alpha1)
  }
  
  model1[i,1]<-bmp1
  model2[i,1]<-bmp2
  model3[i,1]<-bmp3
}


model1<-data.frame(model1)
model2<-data.frame(model2)
model3<-data.frame(model3)


colnames(model1)<-c('Black market premium')
model1$month<-seq(1,dim(model1)[1])
model1<-melt(model1,id.vars = 'month')

colnames(model2)<-c('Black market premium')
model2$month<-seq(1,dim(model2)[1])
model2<-melt(model2,id.vars = 'month')

colnames(model3)<-c('Black market premium')
model3$month<-seq(1,dim(model3)[1])
model3<-melt(model3,id.vars = 'month')

model1$model<-'Model 1'
model2$model<-'Model 2 - Time trend'
model3$model<-'Model 3 - Proportion of term'

results<-rbind(model1,model2,model3)


# Plot

ggplot(data=subset(results,month<=20 ),
       aes(x=month,y=value,color=model,linetype=model))+
  geom_line()+theme_bw()+theme_classic()+
  ylab('Monthly Effect on Presidential Approval')+xlab('Month in which the effect takes place')+
  theme(legend.position = c(.85,.2))+ labs(color="")+
  scale_x_continuous(breaks=seq(1,dim(subset(results,month<=20))[1],1))+
  scale_y_continuous(breaks=seq(-3,0,.25))+
  scale_colour_manual(name = "Specification",
                      labels = c("Model 1",
                                 "Model 2 - Time trend", 
                                 "Model 3 - Proportion of term"),
                      values = c("black","red" ,"blue")) +   
  scale_linetype_manual(name = "Specification",
                        labels = c("Model 1",
                                   "Model 2 - Time trend", 
                                   "Model 3 - Proportion of term"),
                        values = c(1, 2,
                                   3)) 



